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£>y Abstract. We present an ab initio and non-perturbative time-dependent approach to the 

i-£h problem of double ionization of a general atom driven by intense XUV laser pulses. After 

using a highly flexible B-Spline i?-matrix method to generate field-free Hamiltonian and 
electric dipole matrices, the initial state is propagated in time using an efficient Arnoldi- 
Lanczos scheme. Test calculations for double ionization of He by a single laser pulse yield 
' good agreement with benchmark results obtained with other methods. The method is then 

i applied to two-color pump-probe processes, for which momentum and energy distributions of 

the two outgoing electrons are presented. 

o 

Q\ ' PACS numbers: 42.50.Hz, 32.60.+i, 32.80.Rm, 32.80.Fb 

O 

■ 1. Introduction 

The two-photon double ionization (DI) of the helium atom induced by intense short XUV 
laser pulses has received considerable attention from both theorists and experimentalists alike. 
Instead of listing a large number of references here, we note that much of the recent work was 
quoted in recent papers 03 HI. Even within the past few months, however, several additional 
papers appeared. 

Given the intensities and lengths of the laser pulses involved, the numerical approaches 
used to tackle this problem are all essentially attempts to solve the time-dependent 
Schrodinger equation (TDSE), beginning with a well-defined initial state before the laser 
strikes and then propagating this state in the presence of the laser field by one of a number 
of numerical approaches. Once the laser is switched off, various probabilities and, in some 
cases, generalized cross sections can be extracted. 
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Over the past two years, our group has been working on the development of a general 
ab initio theoretical approach, which is applicable to complex targets beyond (quasi) two- 
electron systems. In two recent papers [fJl H) we outlined how field-free Hamiltonian and 
electric dipole matrices generated with the highly flexible 5-Spline .R-matrix (BSR) [5j suite 
of codes may be combined with an efficient Arnoldi-Lanczos time propagation scheme to 
describe the interaction of short intense laser pulses with a complex atom, leading to multi- 
photon excitation and single ionization. The key points of our method are the following: 
1) We employ the BSR code, which allows for the use of non-orthogonal orbital sets, to 
generate field-free Hamiltonian and dipole matrices. 2) We then set up an efficient Arnoldi- 
Lanczos scheme to propagate the initial state in time. 3) Finally, we extract the information 
by standard projection schemes. 

In the present paper, we report on the extension of this approach and the corresponding 
computer code [|6l to allow for two electrons in the continuum and hence the possibility 
to describe double ionization processes. After outlining the general method, we present 
a test application to the He problem, for which many benchmark results are available for 
comparison. Finally, we use the method to investigate pump-probe processes involving two 
XUV laser pulses whose characteristics, including a time delay, are assumed to be controlled 
separately in the corresponding experiment. 

Unless specified otherwise, atomic units (a.u.) are used throughout this manuscript. 

2. THEORY 

2.1. The TDBSR Approach 

The ingredients of an appropriate theoretical and computational formulation require an 
accurate and efficient generation of the Hamiltonian and the electron— field interaction matrix 
elements, as well as an optimal approach to propagate the TDSE in real time. As mentioned 
above, there have been numerous calculations for two-electron systems such as He and also 
H 2 . While these investigations emphasize the important role of two-electron systems in 
studying electron— electron correlation in the presence of a strong laser field, in its presumably 
purest form, experiments with He atoms are difficult and other noble gases, such as Ne and 
Ar, are often favored by the experimental community. 

Fully ab initio theoretical approaches, which are applicable to complex targets beyond 
(quasi) two-electron systems, are still rare. For (infinitely) long interaction times, the R- 
matrix Floquet ansatz has been highly successful. A critical ingredient of this method is 
the general atomic .R-matrix method developed over many years by Burke and collaborators 
in Belfast. A modification of the method, allowing for relatively long though finite-length 
pulses was described by Plummer and Noble [8]. Recently, a time-dependent formulation [|9l 
of this method was applied to short-pulse laser interactions with Ar ifTOl and a pump-probe 
XUV scenario for Ne [fTTIl . 

Following our recent work on excitation and single-ionization of Ne [3] and Ar BH, we 
now describe what needs to be done to extend our approach to double ionization. We start 
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with the time-dependent Schrodinger equation 

d 

i— \P( ri , r N - 1) = [H (r u r N ) + V(r u r^; *p(n, r w ; t) (1) 

for the iV-electron wavefunction ^(r±, . . . , Vn', t), where H (ri, r^) is the field-free 
Hamiltonian containing the kinetic energy of the iV electrons, their potential energy in the 
field of the nucleus, and their mutual Coulomb repulsion, while 

N 

V(r 1 ,...,r N ;t) = J2E(t)-r l (2) 

represents the interaction of the electrons with the laser field E{t) in the dipole length form. 

The tasks to be carried out in order to computationally solve this equation and to extract 
the physical information of interest are: 

1. Generate a representation of the field-free Hamiltonian and its eigenstates; these include 
the initial bound state, other bound states, autoionizing states, as well as single- 
continuum and double-continuum states to represent electron scattering from the residual 
ion. 

2. Generate the dipole matrices to represent the coupling to the laser field. 

3. Propagate the initial bound state until some time after the laser field is turned off. 

4. Extract the physically relevant information from the final state. 

Of particular interest in the experiments mentioned above are processes, in which 
one, two, or even more electrons undergo significant changes in their quantum state in the 
presence of an atomic core. These include excitation, single and double ionization, ionization 
plus simultaneous excitation, or inner-shell ionization with subsequent rearrangement in the 
hollow ion. The latter processes, in particular, can only be investigated in systems beyond 
the frequently studied two-electron helium atom or H 2 molecule. This generalization to two 
electrons outside a multi-electron core is far from trivial, but the flexibility of the BSR method 
is highly advantageous for the tasks 1 and 2. In contrast to many other approaches lfT2l[T3l[T4l . 
our method is formulated in a sufficiently general way to be applicable to complex atoms, such 
as inert gases other than helium and even open-shell systems with non-vanishing spin and 
orbital angular momenta. In reality, of course, the size of the problems that can be handled is 
also determined by the available computational resources. 

The solution of the TDSE requires an accurate and efficient generation of the 
Hamiltonian and electron— field interaction matrix elements. In order to achieve this goal, 
we approximate the time-dependent wavefunction as 

(ri, r N ; t) « C q (t)$ q ( ri , r N ). (3) 
q 

The $g(ri, r N ) are a set of time-independent iV-electron states formed from appropriately 
symmetrized products of atomic orbitals. They are expanded as 

$,(ri, r N ) = a ijcq®c(xi, • • • , xn-2] f N ^a N ^i; r^cr n)Ri (rzv- i ) Rj ( r iv) ■ (4) 
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Here A is the anti-symmetrization operator, the c (xi, . . . , ^n-i^n-i] ^n^n) are 
channel functions involving the space and spin coordinates (Xj) of N — 2 core electrons 
coupled to the angular (f) and spin (a) coordinates of the two outer electrons, Ri(r) is a 
radial basis function, and the a i; j cq are expansion coefficients. Although resembling a close- 
coupling ansatz with two continuum electrons, the expansion §4§ contains bound states and 
singly ionized states as well. In general, the atomic oritals, Ri(r), are not orthogonal to one 
another or to the orbitals used to describe the atomic core. If orthogonality constraints are 
imposed on these functions, additional terms would need to be added to the expansion to 
relax the constraints. This possibility still exists as an option in our computer code. 

In addition to simplifying the expansion, a significant advantage of not forcing such 
orthogonality conditions is the flexibility gained by being able to tailor the optimization 
procedures to the individual neutral, ionic, and continuum orbitals. In the BSR code, the 
outer orbitals (i.e., the R functions above) are expanded in £>-splines. Factors that depend 
on angular and spin momenta are separated from the radial degrees of freedom through 
the construction of the channel functions. Since many Hamiltonian matrix elements share 
common features, this enables the production of a "formula tape", resulting in an efficient 
procedure to generate the required matrix elements. 

When the expansion © is inserted into the Schrodinger equation, we obtain 

iS^C{t) = [H + E(t)D] ■ C(t), (5) 

where S is the overlap matrix of the basis functions, H and D are matrix representations of 
the field-free Hamiltonian and the dipole coupling matrices, and C(t) is the time-dependent 
coefficient vector in Eq. ©. 

The price to pay for the flexibility in the BSR approach, at least initially, is the 
representation of the field-free Hamiltonian and the dipole matrices in a non-orthogonal basis. 
In our previous work Q, we described two methods of combining the BSR method with a 
highly efficient Arnoldi-Lanczos propagation scheme. In the first one, we diagonalized the 
overlap matrix S and transformed the problem back to an orthogonal basis before applying 
the standard propagation scheme [fT5l[T6l . Alternatively, the propagation can be done more 
directly in the non-orthogonal basis. This only requires a Cholesky decomposition of S, but 
some additional operations at every time step. 

A third approach involves a transformation into the eigenbasis by solving the field-free 
generalized eigenvalue problem first. Details can be found in Guan et al (31 . This method, 
also used by Laulan and Bachau ifTTil . has the major advantage that the transformation makes 
it possible to cut unphysically high eigenvalues and the corresponding eigenvectors from the 
time propagation scheme, thereby making the scheme more stable. At the same time, it 
simplifies the definition of the initial state and the extraction of the physically interesting 
information, and it once again allows for the use of the standard Arnoldi-Lanczos time 
propagation scheme. It was applied in the present work as well. We checked that the results 
were stable against a variation of both the cut-off energy (5.5 a.u. for the results shown below) 
and the time-step (100 steps per optical cycle). 

We first tested our method for the He problem, i.e., two electrons in the presence of a 
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bare He 2+ core, and will show some example results below. With a box size of 60 a.u. and 81 
£>-splines distributed over this range, with a spacing of 1.0 a.u. near the edge of the box and 
smaller steps near the origin, the ranks of the field-free Hamiltonian blocks ranged between 
10,000 and 20,000. We expect these ranks to increase to 50,000—100,000 for a complex target 
like Ne. Fortunately, it has become rather straightforward to diagonalize matrices of that size 
and perform the transformation of the dipole matrices to the eigenbasis on massively parallel 
computing platforms. 



2.2. Extraction of the cross sections in two-photon double ionization process 

Details regarding the extraction of the fully differential (energies and angles resolved) as well 
as the total cross section were given by Guan et al and therefore will not be repeated here. 
To make this paper self-consistent, however, we summarize the most important formulas used 
for single-color processes and then indicate what has to be done in the case of more than one 
laser. 

In a time-dependent formulation, the dependence of the iV-photon cross section on the 
number of photons absorbed not only occurs through the factor (u/Io) N , where to is the 
angular frequency and I is the peak intensity of the laser, but also through the effective 
interaction time. The result for the latter depends explicitly on the number of photons being 
absorbed. For two photons it is given by [18] 

t c ( ; } = [ T f 2N (t)dt ( 6) 

Jo 

with the result 

T® = 35r/128 (7) 

for a sine-squared pulse shape of duration of r. 

The energy sharing between the two escaping electrons can be uniquely determined 
through the hyperangle a = tan -1 ^//^). We set E\ = E exc cos 2 a and E 2 = E cxc sin 2 a, 
where E exc = 2u — I 2+ is the excess energy for two-photon double ionization. For a given 
energy sharing a, the triple-differential cross section (TDCS) can be written as 

where \^f(t)) is the time-propagated wavefunction. Also, the singlet two-electron continuum 
wavefunction satisfying the incoming-wave boundary condition (— ) is given by 

1 



<Vi)4;V 2 )+ <Vi)4;V 2 ) 



(9) 



V2 

With the one-electron Coulomb function given by 

Ira 

and the asymptotic behavior 

(c) / \ / 2 , Z In 

Vkiv) — > \/ - sin (fcr + — ln2fcr - — + a;), (11) 
r — > oo V 7T K 2 



Double Ionization of Atoms by XUV Laser Pulses 



6 



with o\ as the Coulomb phase, the two-electron Coulomb function © is normalized in 
momentum space according to 

(*tij*£*i> = ^ " ^l)^ " k 2) + ~ ~ K)- (12) 

In the present work, we first generated the one-electron Coulomb functions <^if ( r ) by 
employing the routine COULFG of Barnett [[T9l , on the same grid used for the 5-splines. This 
choice guarantees a self-consistent grid representation of the time-evolved wavepacket and 
the final continuum states. We then expressed the function in terms of our £>-spline basis in 
order to use the standard integration schemes in that basis. 

It is important to remember that there are two indistinguishable electrons in the final 
channel, and hence the energy sharings described by a and n/2 — a represent the same 
observable event. Therefore, we either need to consider < a < n/Aor n /A < a < 7r/2to 
avoid double counting. The TDCS with respect to the energy of one electron is then given by 
d 3 cx 1 /lo\ 2 I 

d£idfeidfc 2 hk 2 cos 2 a \I J t$ 



dk' 1 dk 2 k[8(k 2 — A^tana) 



x 



iLMi 

L=0,2,ZiZ 2 



(13) 



Here we have defined the two-electron projection 

T^M^^l^m)), d4) 

with the phase factor 

XwhikiM) = (-i^+VK^+^to)), (15) 

and 

<tU = ^{rilCri)^!^)^^!,^)} (16) 
as the properly antisymmetrized product of two one-electron radial Coulomb functions and 
coupled spherical harmonics yh, for the angular part. [Since M = for linearly polarized 
radiation, we omit it to simplify the notation.] 

By collecting all ionization events, the total cross section for two-photon double 
ionization is obtained as 

[a (a* At dV f u Y 1 fAi,Ai,^hM n-n 
ctns = da dfcidfc 2 =1 — 1 _ - / dhdk 2 , (17) 

J J dadfcxdfca VV T£ J dk x dk 2 

where we have introduced the momentum distribution 



d 2 V(h,k 2 ) _ V V| L , 



dk\dk 2 

L=0,2 hl 2 



Similarly, the momentum distribution is given by 

h, k 2 ) 



d 2 S(E 1 ,E 2 ) 1 d 2 V(h,k 2 



(18) 



(19) 



dE 1 dE 2 k\k 2 dk\dk 2 
Below we will also show results obtained for two linearly polarized XUX laser pulses, 
for which the electric field is given by 

E(t) = sin(^) + E 2 f 2 (t - T d ) sm{u 2 {t - T d ))]e z . (20) 
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Here uj\ and lu 2 are the two central frequencies while T d is the time delay between the two 
pulses. The total time duration can be written as r = max{r 1; r 2 + T d }. 

Note that this time duration depends on the durations of pulse 1, pulse 2, and the time 
delay between them. Hence, the definition of an "effective interaction time" is by no means 
straightforward and even a "generalized cross section" (which should be independent of these 
parameters) cannot be defined in the standard sense. Momentum distributions, however, can 
be obtained by calculating and then squaring the magnitude of the function (fT4~l) . Energy 
distributions for the He case were recently reported by Foumouo et al ll2~0l . 

3. RESULTS AND DISCUSSION 

As a first test of our newly developed approach, we applied it to the single-color two-photon 
double ionization of helium, for which a variety of benchmark data exist for comparison (see, 
for example, refs. HI 13 and references therein). In this particular test, we confined the system 
to a spatial box of r max = 60 a.u. We used 81 £>-splines of order 8 to span this configuration 
space, with the knot sequence chosen in such a way that the spatial variation of the wave- 
function close to the nucleus (determined by the nuclear charge Z) and far away from the 
nucleus (determined by the highest energy of a free electron) could be represented accurately. 

Having obtained the initial state of the system by diagonalizing the field-free Hamiltonian 
of the X S C symmetry and taking the eigenvector to the lowest eigenvalue (—2.90330 a.u. in our 
particular case), the time propagation in the laser field is accomplished through the Arnoldi- 
Lanczos algorithm. Compared to other time-propagation approaches, such as leapfrog or 
a split-operator approach [12111 . the present scheme allows us to take relatively large steps 
in time. Specifically, using only 100 steps per optical cycle (o.c.) is sufficient to achieve 
converged solutions of the TDSE for the cases presented in this paper. 

The results presented below were obtained with a sine- squared pulse with a peak intensity 
of 5 x 10 14 W/m 2 and a time duration of 10 optical cycles. For this case, we first checked the 
total cross sections for photon energies of 42 and 50 eV, respectively. These values were 
chosen to ensure that the box size and the length of the laser pulse were suitable to properly 
define a generalized cross section for two-photon double ionization. The present results of 
4.02 x 10~ 53 cm 4 s and 2.20 x 10~ 52 cm 4 s agree in a satisfactory way with the corresponding 
values of 3.76 x 10~ 53 cm 4 s and 1.89 x 10~ 52 cm 4 s obtained in an entirely different finite- 
element discrete- variable representation (FE-DVR) 0. 

3.1. Triple-differential cross sections for one-color two-photon double ionization 

Figure CD displays the triple-differential cross section for two-photon double ionization of 
helium in a 10-cycle sine-squared laser pulse of central photon energy 42 eV and peak 
intensity 5 x 10 14 W/cm 2 for equal energy sharing (E\ = E 2 = 2.5 eV) of the two outgoing 
electrons. These results are for the coplanar geometry, where the electric field vector of the 
linearly polarized laser field and the momentum vectors of the two escaping electrons all lie 
in the same plane. 




Figure 1. The coplanar triple-differential cross section for two-photon double ionization of 
helium in a 10-cycle sine-squared laser pulse of central photon energy 42 eV and peak intensity 
5 x 10 14 W/cm 2 for equal energy sharing (Ei = Ei = 2.5 eV) of the two outgoing electrons. 
The angle listed in the figure is the angle between the laser polarization vector and one of the 
two escaping electrons, while the emission angle of the second electron varies. The present 
TDBSR results are compared with those from our previous FE-DVR calculation Q . 



Given the very different implementations of the current TDBSR method and the FE- 
DVR approach, we see satisfactory agreement between the present results and those obtained 
earlier. Note that these results were obtained with (£ max , / liiriax , Z 2imax ) — (3,3,3). We 
are aware of potential convergence problems with these parameters, especially for kinematic 
situations where the cross section is small (e.g., at 6\ = 90°). However, since the principal 
motivation for performing these particular calculations was to test the method and the 
accompanying computer code, we consider the test successful and hence move on to the 
two-color-case in the next subsection. 

3.2. Momentum and energy distributions for two-color two-photon double ionization 

We now consider the process of double photoionization by absorption of two photons at 
different central photon energies. In other words, the target helium atom is exposed to the 
irradiation by two laser pulses, of potentially different frequencies and with a controllable 
time delay. This allows us to study the mechanism of the breakup problem by applying pulses 
of various durations for each pulse and also modifying the time delay between them. The 
following two-color laser parameters were considered in this work: pulse 1 has a central 
photon energy of 35.3 eV, a time duration of 12 o.c, and a peak intensity of 1 x 10 14 W/cm 2 , 
while the corresponding parameters for pulse 2 are 57.1 eV, 14 o.c, and 1 x 10 13 W/cm 2 , 
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Figure 2. Electric field of the two-color laser pulse. The laser parameters are: lo\ — 35.3 eV, 
Ti = 12 o.c. at a peak intensity of 10 14 W/cm 2 (solid curve), and 0J2 = 57.1 eV, T2 = 14 o.c. at 
a peak intensity of 10 13 W/cm 2 (dashed curve). The time delays between the maxima of the 
two pulses are, from left to right, —16.5, 0.0, and 16.5 a.u. 



respectively. A similar two-color laser pulse was studied in the recent work by Foumouo et al 

US). 

In the present work, we are interested in the mechanism for the ejection of two electrons 
when the time delay between the two pulses is varied. Figure [2] displays the time-dependent 
electric fields of the two-color pulses for time delays of —16.5 a.u. (i.e., pulse 2 comes first), 
no delay (both pulses come simultaneously, and a delay of +16.5 a.u. (pulse 1 comes first). 
Note that +16.5 a.u correspond to approximately 400 attoseconds. 

Figure [3] depicts the momentum and energy distributions of the two escaping electrons 
for the above two-color laser pulses. For the negative time delay with pulse 2 arriving first, 
we observe a strong peak around E\ = E 2 — 6.7 eV, as well as two weaker peaks around 
{Ei, E 2 ) ~ (2.7, 32) eV and (Ex, E 2 ) ~ (32, 2.7) eV, respectively. This suggest that the target 
electrons predominantly absorb the two photons in the highly correlated way similar to the 
nonsequential double ionization mechanism, even though the photon energies are different. 
Recent calculations by Feist et al 11221 . using a single-color laser, also demonstrate that the 
conventional scenario of "direct" vs "sequential" two-photon double ionization breaks down 
for attosecond XUV pulses, even when the sequential process is energetically allowed. 

In the above two-color case, the two ejected electrons share the excess energy of 
E exc = loi + uj 2 — I 2+ — 13.5 eV almost equally. However, there is apparently another 
ionization channel open for the present laser parameters. Since the target electrons interact 
with pulse 2 of energy 57.1 eV first, the first electron is ejected with a kinetic energy of 
about 32 (~ 57.1 — 24.6) eV before the second electron is ionized with an energy of 
2.7 (~ 57.1 — 54.4) eV from the ground state of the residual He + (ls) ion. This channel 
is characteristic of sequential ionization by two photons of the same high frequency and very 
different from the direct channel. In fact, in this case there is no electron left for the second 
pulse to interact with. 

To improve the visibility of the remaining figures, we will show the energy distributions 
only in the energy region below 20 eV for the cases of zero and positive time delays. When 
there is no time delay between the two maxima of the electric fields, the principal structure 
observed is still around E x = E 2 ~ 6.7 eV in the energy distribution. The mechanism in this 
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Figure 3. Momentum (left column) and energy (right column) distributions of the two 
escaping electrons in a two-color laser pulse. The laser parameters are the same as in figure[2] 
The color bars for the first row (negative time delay) correspond to units of 10~ 5 a.u., while the 
other two rows are given in units of 10~ 4 a.u. The time delays between the two laser pulses, 
from top to bottom, are —16.5, 0.0, and +16.5 a.u., respectively. 



case is thus similar to what we observed for the negative time delay. 

On the other hand, when the time delay is increased further, for example up to +16.5 a.u., 
the two electrons are ejected a very different way. In this case, the main structure revealed 
from figure |3] exhibits two peaks at (E u E 2 ) ~ (10.0, 2.5) eV and (2.5, 10.0) eV. This is due 
to the fact that pulse 1 comes first and hence, after absorption of one photon of frequency 
one electron with the kinetic energy of 10.7 (= 35.3 — 24.6) eV is ejected. The time delay 
of +16.5 a.u. is sufficiently long for the residual He + ion to relax into its ground state after 
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Figure 4. Momentum (left) and energy (right) distributions of the two escaping electrons in a 
two-color laser pulse. The wavelengths and laser intensities are the same as in figure[2] but in 
this case the pulselengths are only 10 optical cycles with a delay of 15 a.u. The color bars are 
given in units of 10~ 4 a.u. 



the single ionization. When this process is followed by the absorption of the second photon 
of 57.1 eV the other electron carries away the excess energy of 2.7 = (57.1 — 54.4) eV. 
This suggests that the two electrons are basically removed through a sequential process. 

We then explored the parameter space in more detail by varying both the pulselength and 
the delay. An even more striking example of the competing processes is exhibited in figure |U 
For the same wavelengths and intensities as before, but with a relatively small change in the 
pulselengths to 10 optical cycles and a delay of 15 a.u. 360 attoseconds), sequential and 
non- sequential absorption of the two photons play an almost equally important role in the 
double ionization process. 

From the discussion and examples given above, it is clear that the time delay plays a 
decisive role in determining how the two electrons are ejected by two-color XUV laser pulses. 
Depending on the details of time delay, the electrons can be ejected in ways either similar to 
the sequential or the non- sequential process. Our findings qualitatively agree with those of 
Foumouo et al Il2~0l for the two-color problem and Feist et al [|22| in the single-color problem. 
They serve as an independent confirmation of their predictions, and also provide us confidence 
in our newly developed computer code. 

4. CONCLUSIONS 

We have presented a general method to calculate two-photon double ionization of atoms 
using a 5-spline i?-matrix approach in connection with an efficient Arnoldi-Lanczos time 
propagation scheme. Test calculations for helium revealed good agreement with previous 
benchmark results obtained with different and entirely independent methods. The method was 
then applied to two-photon double ionization of helium using two time-delayed XUV laser 
pulses. The latter results confirm that attosecond spectroscopy will provide a "microscope" to 
examine and also control the way electrons interact in atomic and molecular targets. 

We are currently in the process of generating and then transforming the corresponding 
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matrices for the two-photon double ionization problem of neon atoms. This will allow us to 
make a direct comparison with the recent experiments of Moshammer et al 11231 carried out 
at the FLASH facility in Hamburg. While the additional complication of a residual core with 
non-zero angular momentum is substantial, the current method has been formulated in such 
a way that these calculations are effectively limited by the available hardware (i.e. super- 
computer facilities) rather than special-purpose software. As long as the field-free and 
the dipole matrices can first be generated, in any primitive basis, and then transformed 
to the eigenbasis, the time propagation essentially becomes a matter of fast matrix-vector 
multiplications, for which highly tuned subroutines are available. 

Consequently, we are confident that we will be able to generate results for complex 
targets in the very near future. In addition, we plan to adapt the method to simple linear 
molecules that are amenable to single-center expansion techniques. 
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